\chapter{引力计算}
\section{论均匀球体内质点的引力}
\subsection{李国斌注记2020.03.09}
牛顿在《自然哲学之数学原理》第一编第12章球体的吸引力证明了万有引力定律，即命题76定理36推论3(参见\ref{NewtonsLawofUniversalGravitationMass})推论4(参见\ref{NewtonsLawofUniversalGravitationDistance})。牛顿通过证明定理70、71，得到球壳内质点不受球壳引力作用，球壳外质点受到引力反比于质点到球壳球中心距离平方。牛顿采用均匀球壳积分证明了这些定理，球壳积分的优点是球壳面密度一定是相等的，因为球壳到球心距离相等，尽管球壳与球壳之间密度不等。参看\ref{NewtonGravityOnSphyere}。

但是，下面引用的网络论文作者使用薄板积分，得出了与牛顿明显不同的结论，究竟是作者正确还是牛顿正确？李国斌逐步分析并推导了薄板积分，发现积分似乎没有错误。那么问题到底出在哪里？可能在于\textbf{薄板积分本身。薄板到球心的密度显然是变化的，不能作为常数处理}，这就是错误之所在。而牛顿的球壳积分不存在这个问题，因为球壳到球心距离处处相等，薄球壳密度必然是常数。

按照同样的理由，可以得到结论：静电学中所有假定薄板电荷密度恒定，并以此积分推导与点电荷之间电势、电场强度得到的结果都是违反实际的。

链接：
https://www.docin.com/p-1637042582.html

详细推导见下面链接
https://www.docin.com/p-1474451629.html

经典物理中，牛顿引力公式$F=G\frac{mM}{R^2}$的困境就是如何计算物体内部质点所受力的大小，一般人们用积分求球壳对内部质点引力，认为球壳对内部质点引力为零，以此来断言，内部质点所受力仅仅是质点所在位置到球心部分的球体对于质点的引力，大小为

\begin{equation}
	F(r)=G\frac{mM}{r^2}
\end{equation}
然而笔者计算表明并非如此。以下是计算的结果。
\begin{equation}
	\frac{F}{m}=g(a)=2\pi \rho RG(1+2a-\frac{1-a^2+\frac{a^4}{3}}{(1-a^2)^{\frac{3}{2}}+1})
\end{equation}
或者：
\begin{equation}
	\frac{F}{m}=g(a)=\frac{3MG}{2R^2}(1+2a-\frac{1-a^2+\frac{a^4}{3}}{(1-a^2)^{\frac{3}{2}}+1})
\end{equation}

$a=sin\delta,Rsin\delta$为球体内质点到球心的距离，G为引力常数。

推导过程如下：
\subsection{半径为r的圆形薄片对于圆心轴线上距离为a的质点的引力F}
\begin{figure}
	\caption{圆形薄片引力\label{circlethinpinGravity}}
	\includegraphics[scale=0.3]{circlethinpinGravity}
\end{figure}
设$\rho$为薄片的密度，$m_0$为质点质量，l是质点到薄片的距离,$\Delta d$为薄片的厚度，把圆片沿着径向分成n个均匀宽度为$\Delta r$的同心环，第j个环的质量为$2\pi \rho(j/n)r\Delta r\Delta d$，它对质点$m_0$的引力指向圆环中心，大小为

\begin{align}
	\Delta F&=G m_0 \int_{0}^{2\pi}\rho (j/n)rd\theta\Delta r\Delta d \frac{1}{l^2+((j/n)r)^2}\frac{l}{(l^2+((j/n)r)^2)^{1/2}}\\
	&=Gm_0 2\pi \rho (j/n)r\Delta r\Delta d \frac{l}{(l^2+((j/n)r)^2)^{3/2}}
\end{align}

简记$C=G m_0 2\pi \rho$，那么整个薄片对于质点$m_0$的引力为
\begin{align}
	F&=\lim_{n \to \infty}\sum_{j=0}^{n-1}\Delta F\\
	&=C \lim_{n \to \infty}\sum_{j=0}^{n-1}(j/n)r\Delta r\Delta d \frac{l}{(l^2+((j/n)r)^2)^{3/2}}\label{thinpiegravityEQ1}
\end{align}
令
$z=(j/n)r$，

则
\begin{align}
	dz &= j/n dr \label{thinpiegravityEQ2}\\
	r &=n/j z\label{thinpiegravityEQ3}
\end{align}

代入式\ref{thinpiegravityEQ1}得

\begin{align}
	F&=C \lim_{n \to \infty}\sum_{j=0}^{n-1}z\Delta r\Delta d \frac{l}{(l^2+z^2)^{3/2}} \\
	&=C  \int_{0}^{r}zdz\Delta d \frac{l}{(l^2+z^2)^{3/2}} \\
	&=C \Delta d l \int_{0}^{r} \frac{zdz}{(l^2+z^2)^{3/2}} \\
	&=C \Delta d l  \frac{-1}{(l^2+z^2)^{1/2}} |_{0}^{r}\\
	&=C  (n/j)\Delta d l[1/l-\frac{1}{(l^2+((j/n)r)^2)^{1/2}}]\label{thinpiegravityEQ4}\\
	&=C  \Delta d l(1-\frac{1}{\sqrt{l^2+r^2}}) \label{thinpiegravityEQ5}\\
	&=C  \Delta d (1-\frac{l}{\sqrt{l^2+r^2}}) \label{thinpiegravityEQ6}
\end{align}

即
\begin{equation}
	F=2\pi m_0 \rho G\Delta dl(1-\frac{1}{\sqrt{l^2+r^2}} )\label{thinpiegravityEQ}
\end{equation}

检查上面的推导，发现式\ref{thinpiegravityEQ5}括号后第一项1是错的，应为1/l，因此式\ref{thinpiegravityEQ}也是错的。改正后，式\ref{thinpiegravityEQ6}和下式\ref{thinpiegravityEQ7}才是正确的：

\begin{equation}
	F=2\pi m_0 \rho G\Delta d(1-\frac{l}{\sqrt{l^2+r^2}} )\label{thinpiegravityEQ7}
\end{equation}

\subsection{求均匀球体内部直径方向上某质点受到的引力}
\begin{figure}
	\caption{均匀球体内部径向上某质点受到的引力\label{GravitationOfParticleInSphere}}
	\includegraphics[scale=0.6]{GravitationOfParticleInSphere}
\end{figure}
根据以上原理我们可以把球体分成厚度为$\Delta d$的2k个均匀圆形薄片,假设球半径为R，质点位于x轴，那么有
\begin{align}
	x^2+y^2+z^2&=R^2\\
	\Delta d k&=R\\
	\Delta d&=\Delta x
\end{align}

代入式\ref{thinpiegravityEQ}，得

$\Delta F(x)=2\pi \rho G\Delta x(1-\frac{l}{\sqrt{l^2+r^2}} )$

球体内部点a处受大小半球引力作用，该引力等于两个矢量力F1、F2的和。

而此处对于F1，F2来说l=x-b，r=y
\begin{align}
	F_1&=C \int_{b}^{R}(1-\frac{x-b}{\sqrt{(x-b)^2+y^2}} )dx\\
	F_2&=C \int_{-R}^{b}(1-\frac{x-b}{\sqrt{(x-b)^2+y^2}} )dx
\end{align}

设
\begin{align}
	x=Rcos\alpha\\
	b=Rcos\delta\\
	y=Rsin\alpha
\end{align}


$F=F_2-F_1$
设
$b=\frac{a^2+1}{2a},a=sin\delta$
现在以$sin\alpha$自变量改写积分

\begin{align}
	\mbox{令}cos(\alpha+\phi)&=1+cos^2\delta-2cos\delta cos\alpha\\
	\mbox{则}\cos^2\alpha-sin^2\phi&=1+cos^2\delta-2cos\delta cos\alpha\\
	\cos^2\alpha+2cos\delta cos\alpha-sin^2\phi-1-cos^2\delta&=0\\
	cos\alpha&=\frac{-2cos\delta\pm2\sqrt{2cos^2\delta+1+sin^2\phi}}{2}\\
	cos\alpha&=-cos\delta\pm\sqrt{2cos^2\delta+1+sin^2\phi}\label{solutionTemp1}\\
	sin^2\phi&=\cos^2\alpha+2cos\delta cos\alpha-1-cos^2\delta\\
	sin\phi&=\sqrt{\cos^2\alpha+2cos\delta cos\alpha-1-cos^2\delta}\label{solutionTemp2}
\end{align}

\begin{align}
	F_1&=C \int_{\delta}^{0}(R(cos\alpha-cos\delta)-\frac{(cos\alpha-cos\delta)}{\sqrt{(cos\alpha-cos\delta)^2+sin^2\alpha}} )d(Rcos\alpha)\\
	&=C (\frac{1}{2}R^2(cos\alpha-cos\delta)^2|_{\delta}^{0}-\int_{\delta}^{0}\frac{(cos\alpha-cos\delta)}{\sqrt{(cos\alpha-cos\delta)^2+sin^2\alpha}}d(Rcos\alpha))\\
	&=C (\frac{1}{2}R^2(1-cos\delta)^2-R\int_{\delta}^{0}(\frac{cos\alpha}{\sqrt{1+cos^2\delta-2cos\delta cos\alpha}}dcos\alpha-\frac{cos\delta}{\sqrt{1+cos^2\delta-2cos\delta cos\alpha}}dcos\alpha))\\
	F_2&=C \int_{\pi}^{\delta}(R(cos\alpha-cos\delta)-\frac{(cos\alpha-cos\delta)}{\sqrt{(cos\alpha-cos\delta)^2+sin^2\alpha}} )d(Rcos\alpha)\\
	F_2+F_1&=C(\frac{1}{2}R^2(cos\alpha-cos\delta)^2|_{\pi}^{0}-R\int_{\pi}^{0}(\frac{cos\alpha}{\sqrt{1+cos^2\delta-2cos\delta cos\alpha}}dcos\alpha-\frac{cos\delta}{\sqrt{1+cos^2\delta-2cos\delta cos\alpha}}dcos\alpha))\\
	&=C (\frac{1}{2}R^2((1-cos\delta)^2-(-1-cos\delta)^2)-R\int_{\pi}^{0}(\frac{cos\alpha}{\sqrt{1+cos^2\delta-2cos\delta cos\alpha}}dcos\alpha-\frac{cos\delta}{\sqrt{1+cos^2\delta-2cos\delta cos\alpha}}dcos\alpha))\\
	&=C (\frac{1}{2}R^2(-4cos\delta)-R\int_{\pi}^{0}(\frac{cos\alpha}{\sqrt{1+cos^2\delta-2cos\delta cos\alpha}}dcos\alpha-\frac{cos\delta}{\sqrt{1+cos^2\delta-2cos\delta cos\alpha}}dcos\alpha))\\
	F_1&=C (1-a-\int_{a}^{1}\frac{x-a}{\sqrt{(a^2+1-2ax)^2}} )dx\\
	F_2&=C (1+a-\int_{-1}^{a}\frac{x-a}{\sqrt{(a^2+1-2ax)^2}} )dx\\
	F_1&=C R(1-a-\frac{1}{\sqrt{2a}}\int_{a}^{1}\frac{x-a}{\sqrt{b-x}} )dx\\
	F_2&=C R(1+a-\frac{1}{\sqrt{2a}}\int_{-1}^{a}\frac{x-a}{\sqrt{b-x}} )dx
\end{align}

令
$u=\sqrt{b-x}$，则
\begin{align}
	F_1&=C R(1-a-\frac{1}{\sqrt{2a}}\int_{a}^{1}\frac{b-u^2-a}{u} )d(b-u^2)\\
	F_1&=C R(1-a+\frac{1}{\sqrt{2a}}\int_{a}^{1}b-a-u^2du\\
	F_2&=C R(1+a+\frac{1}{\sqrt{2a}}\int_{-1}^{a}b-a-u^2du\\
	F_1&=C R(1-a+\frac{1}{\sqrt{2a}}[(b-a)u-u^3/3]_{\sqrt{b-a}}^{\sqrt{b-1}}\label{semisphereGravity1}\\
	F_2&=C R(1+a+\frac{1}{\sqrt{2a}}[(b-a)u-u^3/3]_{\sqrt{b+1}}^{\sqrt{b-a}}\label{semisphereGravity1}\\
	F&=C R(2a+\frac{1}{\sqrt{2a}}([(b-a)u-u^3/3]_{\sqrt{b+1}}^{\sqrt{b-a}}+[(b-a)u-u^3/3]^{\sqrt{b-a}}_{\sqrt{b-1}})\\
	F&=C R(2a+\frac{4(b-a)^{\frac{3}{2}}}{3(2a)^2}-\frac{b-a}{a}+\frac{2+6a^2}{3(2a)^2})\\
	F&=C R(2a+\frac{4(1-a)^{\frac{3}{2}}}{3\sqrt{2a}}-\frac{1}{\sqrt{2a}}([(b-a)((b+1)^{1/2}+(b-1)^{1/2})-((b+1)^{3/2}+(b-1)^{3/2})/3])\\	
	F&=C R(1+2a+\frac{(1-a^2)^{3/2}-1}{3a^2})\\
	F&=C R(1+2a-\frac{1-a^2+a^4/3}{(1-a^2)^{3/2}+1})
\end{align}
即
\begin{equation}
	\frac{F}{m_0}=2\pi \rho R G(1+2a-\frac{1-a^2+a^4/3}{(1-a^2)^{3/2}+1})
\end{equation}

$a=sin\delta,Rsin\delta$是质点到球心的距离。
设
$M=\frac{4}{3}\pi \rho R^3$
\begin{equation}
	g(a)=2\pi \rho R G(1+2a-\frac{1-a^2+a^4/3}{(1-a^2)^{3/2}+1})
\end{equation}
上式可改写为
\begin{equation}
	g(a)=\frac{3MG}{2R^2}(1+2a-\frac{1-a^2+a^4/3}{(1-a^2)^{3/2}+1})
\end{equation}

上面推导过程中有存在

\begin{align}
	F_1+F_2&=m_0\frac{16}{3}\pi \rho RG\\
	F_2-F_1&=m_0g(a)
\end{align}

球体内引力存在趋近球体表面内侧它所受到的引力最大，中心附近最小接近最大值的
3/16倍，这和人们预言的与距离成正比的关系有所不同，不过它很接近弹性的线性特	征，当然如果正处球心质点各部分受到强度相等方向不同的各向引力，最终质点受力为零。球心为平衡点，物质无法达到这一点，如果能达到那么它的体积必须为零。计算结果可以发现均匀球体球内质点与球外质点
引力存在很大差异，而且球面是不连续点，在趋近球外壳内侧的引力是趋近球壳外侧的引力的4倍，而且球心区域是外壳外侧的引力0.75倍，而并非我们想象的零...

在作者的另外一篇类似文章中，给出了另外一个计算结果。作者写道：

球体内部引力计算结果为：
\begin{align}
	F&=\frac{4}{3}\pi \rho G(3\sqrt{1-a^2-\frac{3-4a^4}{(1+2a^2)\sqrt{1-a^2}+1}})\\
	g(a)&=g_0(3\sqrt{1-a^2-\frac{3-4a^4}{(1+2a^2)\sqrt{1-a^2}+1}})\\	
	g_0&=\frac{4}{3}\pi \rho G
\end{align}

引力函数是个偶函数，内部由引力产生的压力也可以计算其值为F2-F1，此处未计算，压力接近线性。如果是地球就需要考虑实际相对密度的变化，积分时密度函数就不能是常数，以下是按照地球模型在把密度按不同密度与半径的球体叠加后的模型引力图。

\begin{figure}
	\includegraphics[scale=0.6]{earthgravityModel}
	\caption{地球引力模型图\label{earthgravityModel}}
\end{figure}

核心区引力非零，非但如此，核心是地球引力最强的区域。·它的平衡由物质的核，电等力去平衡。

一般往往会把球面积分当成径向引力积分，错误就会由此产生...


外部引力计算结果为
\begin{equation}
	F=2\pi\rho G\int_{1}^{-1}(1-\frac{a-x}{\sqrt{r^2+(a-x)^2}})dx=2\pi\rho G\int_{1}^{-1}(1-\frac{a-x}{\sqrt{1+a^2-2ax}})dx\\
	\frac{F}{2\pi\rho G}=
\end{equation}

(ligb 2019-06-10感觉作者推导过程有很多错误，不敢抄下去了。)

\subsection{求均匀半球体对大圆面边界直径中心质点的引力}
根据公式\ref{semisphereGravity1}，得
\begin{align}
	F_1&=C R(1-a+\frac{1}{\sqrt{2a}}[(b-a)u-u^3/3]_{\sqrt{b-a}}^{\sqrt{b-1}}\label{semisphereGravity3}\\
	F_1&=C R(1-a+\frac{1}{\sqrt{2a}}[(b-1)^{1/2}(b-a-(b-1)/3)-(b-a)^{1/2}(b-a-(b-a)/3)]\\
	F_1&=C R(1-a+\frac{1}{\sqrt{2a}}[(b-1)^{1/2}((2b-3a+1)/3)-(b-a)^{1/2}(2(b-a)/3)]\\
	F_1&=C R(1-a+\frac{1}{\sqrt{2a}}[(b-1)^{1/2}((2b-3a+1)/3)-2/3(b-a)^{3/2}]
\end{align}
\subsection{求均匀球体内部直径方向上某质点受到的引力}

\begin{figure}
	\caption{均匀球体内部径向上某质点受到的引力\label{GravitationOfParticleInSphere}}
	\includegraphics[scale=0.6]{GravitationOfParticleInSphere}
\end{figure}
根据以上原理我们可以把球体分成厚度为$\Delta d$的2k个均匀圆形薄片,假设球半径为R，质点位于x轴，那么有
\begin{align}
	x^2+y^2+z^2&=R^2\\
	\Delta d k&=R\\
	\Delta d&=\Delta x
\end{align}

代入式\ref{thinpiegravityEQ}，得

$\Delta F(x)=2\pi \rho G\Delta x(l-\frac{l}{\sqrt{l^2+r^2}} )$

球体内部点a处受大小半球引力作用，该引力等于两个矢量力F1、F2的差。

而此处对于F1，F2来说l=x-b，r=y
\begin{align}
	F_1&=C \int_{b}^{R}(x-b-\frac{x-b}{\sqrt{(x-b)^2+y^2}} )dx\\
	F_2&=C \int_{-R}^{b}(x-b-\frac{x-b}{\sqrt{(x-b)^2+y^2}} )dx
\end{align}

设$x=Rsin\alpha,b=Rsin\delta$

$F=F_2-F_1$
设
$b=\frac{a^2+1}{2a},a=sin\delta$
现在以$sin\alpha$自变量改写积分
\begin{align}
	F_1&=C R\int_{a}^{1}(1-\frac{x-a}{\sqrt{(a^2+1-2ax)^2}} )dx\\
	F_2&=C R\int_{-1}^{a}(1-\frac{x-a}{\sqrt{(a^2+1-2ax)^2}} )dx\\
	F_1&=C (1-a-\int_{a}^{1}\frac{x-a}{\sqrt{(a^2+1-2ax)^2}} )dx\\
	F_2&=C (1+a-\int_{-1}^{a}\frac{x-a}{\sqrt{(a^2+1-2ax)^2}} )dx\\
	F_1&=C R(1-a-\frac{1}{\sqrt{2a}}\int_{a}^{1}\frac{x-a}{\sqrt{b-x}} )dx\\
	F_2&=C R(1+a-\frac{1}{\sqrt{2a}}\int_{-1}^{a}\frac{x-a}{\sqrt{b-x}} )dx
\end{align}

令
$u=\sqrt{b-x}$，则
\begin{align}
	F_1&=C R(1-a-\frac{1}{\sqrt{2a}}\int_{a}^{1}\frac{b-u^2-a}{u} )d(b-u^2)\\
	F_1&=C R(1-a+\frac{1}{\sqrt{2a}}\int_{a}^{1}b-a-u^2du\\
	F_2&=C R(1+a+\frac{1}{\sqrt{2a}}\int_{-1}^{a}b-a-u^2du\\
	F_1&=C R(1-a+\frac{1}{\sqrt{2a}}[(b-a)u-u^3/3]_{\sqrt{b-a}}^{\sqrt{b-1}}\label{semisphereGravity1}\\
	F_2&=C R(1+a+\frac{1}{\sqrt{2a}}[(b-a)u-u^3/3]_{\sqrt{b+1}}^{\sqrt{b-a}}\label{semisphereGravity1}\\
	F&=C R(2a+\frac{1}{\sqrt{2a}}([(b-a)u-u^3/3]_{\sqrt{b+1}}^{\sqrt{b-a}}+[(b-a)u-u^3/3]^{\sqrt{b-a}}_{\sqrt{b-1}})\\
	F&=C R(2a+\frac{4(b-a)^{\frac{3}{2}}}{3(2a)^2}-\frac{b-a}{a}+\frac{2+6a^2}{3(2a)^2})\\
	F&=C R(2a+\frac{4(1-a)^{\frac{3}{2}}}{3\sqrt{2a}}-\frac{1}{\sqrt{2a}}([(b-a)((b+1)^{1/2}+(b-1)^{1/2})-((b+1)^{3/2}+(b-1)^{3/2})/3])\\	
	F&=C R(1+2a+\frac{(1-a^2)^{3/2}-1}{3a^2})\\
	F&=C R(1+2a-\frac{1-a^2+a^4/3}{(1-a^2)^{3/2}+1})
\end{align}
即
\begin{equation}
	\frac{F}{m_0}=2\pi \rho R G(1+2a-\frac{1-a^2+a^4/3}{(1-a^2)^{3/2}+1})
\end{equation}

$a=sin\delta,Rsin\delta$是质点到球心的距离。
设
$M=\frac{4}{3}\pi \rho R^3$
\begin{equation}
	g(a)=2\pi \rho R G(1+2a-\frac{1-a^2+a^4/3}{(1-a^2)^{3/2}+1})
\end{equation}
上式可改写为
\begin{equation}
	g(a)=\frac{3MG}{2R^2}(1+2a-\frac{1-a^2+a^4/3}{(1-a^2)^{3/2}+1})
\end{equation}

上面推导过程中有存在

\begin{align}
	F_1+F_2&=m_0\frac{16}{3}\pi \rho RG\\
	F_2-F_1&=m_0g(a)
\end{align}

球体内引力存在趋近球体表面内侧它所受到的引力最大，中心附近最小接近最大值的
3/16倍，这和人们预言的与距离成正比的关系有所不同，不过它很接近弹性的线性特征，当然如果正处球心质点各部分受到强度相等方向不同的各向引力，最终质点受力为零。球心为平衡点，物质无法达到这一点，如果能达到那么它的体积必须为零。计算结果可以发现均匀球体球内质点与球外质点引力存在很大差异，而且球面是不连续点，在趋近球外壳内侧的引力是趋近球壳外侧的引力的4倍，而且球心区域是外壳外侧的引力0.75倍，而并非我们想象的零...

在作者的另外一篇类似文章中，给出了另外一个计算结果。作者写道：

球体内部引力计算结果为：
\begin{align}
	F&=\frac{4}{3}\pi \rho G(3\sqrt{1-a^2-\frac{3-4a^4}{(1+2a^2)\sqrt{1-a^2}+1}})\\
	g(a)&=g_0(3\sqrt{1-a^2-\frac{3-4a^4}{(1+2a^2)\sqrt{1-a^2}+1}})\\	
	g_0&=\frac{4}{3}\pi \rho G
\end{align}

引力函数是个偶函数，内部由引力产生的压力也可以计算其值为F2-F1，此处未计算，压力接近线性。如果是地球就需要考虑实际相对密度的变化，积分时密度函数就不能是常数，以下是按照地球模型在把密度按不同密度与半径的球体叠加后的模型引力图。

\begin{figure}
	\includegraphics[scale=0.6]{earthgravityModel}
	\caption{地球引力模型图\label{earthgravityModel}}
\end{figure}

核心区引力非零，非但如此，核心是地球引力最强的区域。·它的平衡由物质的核，电等力去平衡。

一般往往会把球面积分当成径向引力积分，错误就会由此产生...


外部引力计算结果为
\begin{equation}
	F=2\pi\rho G\int_{1}^{-1}(1-\frac{a-x}{\sqrt{r^2+(a-x)^2}})dx=2\pi\rho G\int_{1}^{-1}(1-\frac{a-x}{\sqrt{1+a^2-2ax}})dx\\
	\frac{F}{2\pi\rho G}=
\end{equation}

(ligb 2019-06-10感觉作者推导过程有很多错误，不敢抄下去了。)
\subsection{质点与均匀球体间的引力}
https://zhidao.baidu.com/question/940618670691950172.html

设质点质量为m，与球心的距离为R。球的半径为a，密度为$\rho$，质量为$M=\rho\frac{4}{3}\pi a^3$。直角坐标系和球坐标系坐标原点都在球心，在两个坐标系质点坐标分别为(0,0,R)和$(R,0,0)$，球内任意微元球坐标为$(r,\phi,\theta)$，其体积为$dV=r^2sin\phi dr d\phi d\theta$，微元和质点的距离为u，建立如图\ref{MassPointGravity}所示的坐标系。

\begin{figure}
	\centering
	\epsfig{file=MassPointGravity,width=8cm}
	\caption{质点与均匀球体间的引力\label{MassPointGravity}}
\end{figure}

根据对称性可知，球对质点的引力必沿z方向，x，y方向上合力为0。

微元对质点的引力为dF

\begin{align}
	C&=Gm\rho\label{MassPointGravityEQ100}\\
	dF&=G\frac{m\rho r^2sin\phi}{u^2}dr d\phi d\theta (R>a)\label{MassPointGravityEQ101}
\end{align}

该引力在z方向上的分力为：

\begin{align}
	dF_z&=dFcos\alpha \label{MassPointGravityEQ102}
\end{align}

由余弦定理，有：

\begin{align}
	u^2&=r^2+R^2-2rRcos\phi\label{MassPointGravityEQ103}\\
	r^2&=u^2+R^2-2uRcos\alpha\label{MassPointGravityEQ104}
\end{align}

式\ref{MassPointGravityEQ103}加式\ref{MassPointGravityEQ104}，得

\begin{align}
	u^2+r^2&=r^2+u^2+2R^2-2rRcos\phi-2uRcos\alpha\\
	0&=R-rcos\phi-ucos\alpha\\
	ucos\alpha&=R-rcos\phi\label{MassPointGravityEQ105}\\
	cos\alpha&=\frac{R-rcos\phi}{u}\label{MassPointGravityEQ106}
\end{align}

对式\ref{MassPointGravityEQ103}两边求微分，得

\begin{align}
	udu&=rRsin\phi d\phi\label{MassPointGravityEQ107}
\end{align}

由式\ref{MassPointGravityEQ103}，得

\begin{align}
	cos\phi&=\frac{r^2+R^2-u^2}{2rR}\label{MassPointGravityEQ108}
\end{align}

式\ref{MassPointGravityEQ106}代入式\ref{MassPointGravityEQ102}，得

\begin{align}
	dF_z&=G\frac{m\rho r^2(R-rcos\phi)sin\phi}{u^3}dr d\phi d\theta\label{MassPointGravityEQ109}\\
	F_s&=\iiint_V dF_z\\
	F_s&=\int_0^{2\pi}\int_0^\pi\int_0^a G\frac{m\rho r^2sin\phi dr}{u^2}\frac{(R-rcos\phi)}{u} d\phi d\theta\label{MassPointGravityEQ110}
\end{align}

式\ref{MassPointGravityEQ107},\ref{MassPointGravityEQ108}代入到\ref{MassPointGravityEQ110}，将自变量$\phi$的函数都替换成u的函数

\begin{align}
	F_s&=Gm\rho\int_0^{2\pi}\int_0^\pi\int_0^a\frac{r^2sin\phi d\phi dr}{u^2}\frac{(R-rcos\phi)}{u} d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^\pi\int_0^a\frac{rRrsin\phi d\phi dr}{Ru^2}\frac{R-rcos\phi}{u} d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\int_{R-r}^{r+R}\frac{rudu dr}{Ru^2}\frac{R-r\frac{r^2+R^2-u^2}{2rR}}{u} d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\int_{R-r}^{r+R}\frac{r}{R}\frac{ududr}{u^2}\frac{2R^2-r^2-R^2+u^2}{2Ru} d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\int_{R-r}^{r+R}\frac{r}{2R^2}\frac{R^2-r^2+u^2}{u^2}dudr d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(\frac{r^2-R^2}{u}+u)|_{R-r}^{r+R}rdr d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(r-R+r+R+(r+R)-(R-r))rdr d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\frac{2}{R^2}r^2dr d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}d\theta\frac{2}{R^2}\frac{1}{3}r^3|_0^a \\
	F_s&=Gm\rho\int_0^{2\pi}\frac{2}{R^2}\frac{1}{3}r^3|_0^a \\
	F_s&=Gm\rho\frac{4\pi}{3}a^3\frac{1}{R^2}
\end{align}
\subsection{质点与均匀球体间的引力解法2}
以下是网络论文提供的另外一种推导方法。
https://wenku.baidu.com/view/e4519371915f804d2b16c1e8.html
\subsection{均匀球体与球体内一个质点间的引力ligb解法}
下面是笔者的推导。

设质点质量为m，与球心A的距离为R。球的半径为a，R<a，密度为$\rho$，球质量为$M=\rho\frac{4}{3}\pi a^3$。质点处B所在球壳包围的球内质量为$M_0=\rho\frac{4}{3}\pi R^3$。

AB直线延长线交球面于南极点O，交球面于北极点N，取O点为直角坐标系和球坐标系坐标原点，ON轴为Z轴且ON方向为Z轴正方向。在两个坐标系中球心坐标分别为(0,0,a)和(a,0,0)，质点坐标为(0,0,a-R),(a-R,0,0)。球面任意微元球坐标为(r,\phi,\theta)，建立如图\ref{MassPointGravity2}所示的坐标系。

质点和球的引力是一种宏观效果，是一个三重积分：

\begin{align}
	C&=Gm\rho\label{MassPointGravityEQ200}\\
	F_s&=Gm\iiint_V\frac{1}{r^2}\rho r^2sin\phi\cos\phi drd\phi d\theta\label{MassPointGravityEQ201}\\
	F_s&=Gm\rho\iiint_Vsin\phi\cos\phi drd\phi d\theta\label{MassPointGravityEQ202}\\
	F_s&=C\iiint_Vsin\phi\cos\phi drd\phi d\theta\label{MassPointGravityEQ203}
\end{align}

剩下就是确定积分限问题。$\theta$上下限显然是从0到$2\pi$。$\phi$区间为[0,$\pi$]。分成两个区间A,B，分别是[0,$\pi$/2],$[\pi/2,\pi]$。在这两个区间，r的上下限不同。

\begin{align}
	\phi_1&=0\label{MassPointGravityEQ204}\\
	\phi_2&=\pi\label{MassPointGravityEQ205}
\end{align}

比较复杂的是确定r的上下限，可以根据几何关系确定：

\begin{align}
	x^2+y^2+(z-R)^2&=(rsin\phi\cos\theta)^2+(rsin\phi\sin\theta)^2+(rcos\phi-R)^2\\
	&=r^2\sin^2\phi+r^2\cos^2\phi-2Rr\cos\phi+R^2=r^2-2Rr\cos\phi+R^2=a^2\\
	r^2-2R\cos\phi r+R^2-a^2&=0\label{MassPointGravityEQ206}
\end{align}

方程\ref{MassPointGravityEQ206}是关于r的一元二次方程，解得

\begin{align}
	r_1&=Rcos\phi-\sqrt{a^2+R^2\sin^2\phi}\label{MassPointGravityEQ207}\\
	r_2&=Rcos\phi+\sqrt{a^2+R^2\sin^2\phi}\label{MassPointGravityEQ208}
\end{align}

两个根$r_1,r_2$划分r的积分区间为$[r_1,r_2]$。

因此，质点和球的引力为

\begin{align}
	F_s&=C\int_0^{\2p}\int_0^{\pi}\int_{r_1}^{r_2}sin\phi\cos\phi drd\phi d\theta\label{MassPointGravityEQ209}
\end{align}

分步积分：

\begin{align}
	\int_0^{\pi}\int_{r_1}^{r_2}sin\phi\cos\phi drd\phi&=\int_0^{\pi}(r_2-r_1)sin\phi\cos\phi d\phi \\
	&=\int_0^{\pi}2\sqrt{a^2+R^2\sin^2\phi}sin\phi\cos\phi d\phi \\
	&=\int_0^{\pi}-\sqrt{a^2+R^2(1+\cos2\phi)/2}d((1+cos2\phi)/2)\\
	&=\int_0^{\pi}-\frac{1}{R^2}\sqrt{a^2+R^2(1+\cos2\phi)/2}d(a^2+R^2(1+cos2\phi)/2)\\
	&=-\frac{1}{R^2}\frac{1}{1/2+1}(a^2+R^2(1+\cos2\phi)/2)^{3/2}|_0^{\pi}\\
	&=-\frac{1}{R^2}\frac{2}{3}(a^2+R^2(1+\cos2\phi)/2)^{3/2}|_0^{\pi}\\
	&=\frac{1}{3}2\sqrt{a^2+R^2(1+\cos2\phi)/2}((2Rcos\phi)^2-R^2cos^2\phi+a^2+R^2\sin^2\phi)\\
	&=\frac{1}{3}2\sqrt{a^2+R^2(1+\cos2\phi)/2}(3R^2(1+cos2\phi)/2+a^2+R^2(1+\cos2\phi)/2)\\
	&\int_0^{\pi}\frac{1}{3}2\sqrt{a^2+R^2(1+\cos2\phi)/2}(3R^2(1+cos2\phi)/2+a^2+R^2(1+\cos2\phi)/2)sin\phi\cos\phi d\phi\\
	&=\int_0^{\pi}\frac{1}{3}\sqrt{a^2+R^2(1+\cos2\phi)/2}(3R^2(1+cos2\phi)/2+a^2+R^2(1+\cos2\phi)/2)sin2\phi d\phi\\
	&=\int_0^{\pi}\frac{1}{3}\sqrt{a^2+R^2(1+\cos2\phi)/2}(3R^2(1+cos2\phi)/2+a^2+R^2(1+\cos2\phi)/2) (-1)d((1+cos2\phi)/2)\\
	&=\frac{-1}{3}(3R^2\frac{2}{15R^4}(3R^2(1+\cos2\phi)/2-2a^2)(a^2+R^2(1+\cos2\phi)/2)^{3/2}+\frac{2}{5R^2}(a^2+R^2(1+\cos2\phi)/2)^{5/2})|_0^\pi=0\\
	F_s&=0
\end{align}

可见均匀球体对球内质点的引力为0。
\subsection{质点与球体间的引力}
设质点质量为m，与球心的距离为R。球的半径为a，球密度为p(r)，质量为$M=p\frac{4}{3}\pi a^3$。直角坐标系和球坐标系坐标原点都在球心，在两个坐标系质点坐标分别为(0,0,R)和$(R,0,0)$，球内任意微元球坐标为$(r,\phi,\theta)$，其体积为$dV=r^2sin\phi dr d\phi d\theta$，微元和质点的距离为u，建立如图\ref{MassPointGravity}所示的坐标系。

\begin{figure}
	\centering
	\epsfig{file=MassPointGravity,width=8cm}
	\caption{质点与球体间的引力\label{MassPointGravity}}
\end{figure}

由于密度p=p(r)仅与r有关，即在一个薄球壳内密度相同，每个薄球壳质量是均匀的。根据对称性可知，球对质点的引力必沿z方向，x，y方向上合力为0。

微元对质点的引力为dF

\begin{align}
	dF&=G\frac{mp(r) r^2sin\phi}{u^2}dr d\phi d\theta (R>a)\label{MassPointGravityEQ101}
\end{align}

该引力在z方向上的分力为：

\begin{align}
	dF_z&=dFcos\alpha \label{MassPointGravityEQ102}
\end{align}

由余弦定理，有：

\begin{align}
	u^2&=r^2+R^2-2rRcos\phi\label{MassPointGravityEQ103}\\
	r^2&=u^2+R^2-2uRcos\alpha\label{MassPointGravityEQ104}
\end{align}

式\ref{MassPointGravityEQ103}加式\ref{MassPointGravityEQ104}，得

\begin{align}
	u^2+r^2&=r^2+u^2+2R^2-2rRcos\phi-2uRcos\alpha\\
	0&=R-rcos\phi-ucos\alpha\\
	ucos\alpha&=R-rcos\phi\label{MassPointGravityEQ105}\\
	cos\alpha&=\frac{R-rcos\phi}{u}\label{MassPointGravityEQ106}
\end{align}

对式\ref{MassPointGravityEQ103}两边求微分，得

\begin{align}
	udu&=rRsin\phi d\phi\label{MassPointGravityEQ107}
\end{align}

由式\ref{MassPointGravityEQ103}，得

\begin{align}
	cos\phi&=\frac{r^2+R^2-u^2}{2rR}\label{MassPointGravityEQ108}
\end{align}

式\ref{MassPointGravityEQ106}代入式\ref{MassPointGravityEQ102}，得

\begin{align}
	dF_z&=G\frac{mp(r) r^2(R-rcos\phi)sin\phi}{u^3}dr d\phi d\theta\label{MassPointGravityEQ109}\\
	F_s&=\iiint_V dF_z\\
	F_s&=\iiint_V G\frac{mp(r) r^2sin\phi dr}{u^2}\frac{(R-rcos\phi)}{u} d\phi d\theta\label{MassPointGravityEQ110}
\end{align}

式\ref{MassPointGravityEQ107},\ref{MassPointGravityEQ108}代入到\ref{MassPointGravityEQ110}，将自变量$\phi$的函数都替换成u的函数

\begin{align}
	F_s&=\iiint_V Gmp(r) \frac{rRrsin\phi d\phi dr}{Ru^2}\frac{R-rcos\phi}{u} d\theta\\
	F_s&=\iiint_V Gmp(r)\frac{rudu dr}{Ru^2}\frac{R-r\frac{r^2+R^2-u^2}{2rR}}{u} d\theta\\
	F_s&=\iiint_V Gmp(r)\frac{r}{R}\frac{ududr}{u^2}\frac{2R^2-r^2-R^2+u^2}{2Ru} d\theta\\
	F_s&=\iiint_V Gmp(r)\frac{r}{2R^2}\frac{R^2-r^2+u^2}{u^2}dudr d\theta\\
	F_s&=\iiint_V Gmp(r)\frac{r}{2R^2}(\frac{R^2-r^2}{u^2}du+du)dr d\theta\\
	F_s&=\iiint_V Gmp(r)\frac{r}{2R^2}(\frac{r^2-R^2}{u}+u)drd\theta\label{MassPointGravityEQ111}
\end{align}

式\ref{MassPointGravityEQ111}是质点和球体之间的通用引力方程。根据该式，可以计算质点与各种球体、非球体的引力。球体密度可以是均匀的(p=常数)，也可以是与半径有关的(p=p(r))。非球体可以是薄或厚球壳，也可以是薄圆环或薄圆板。质点可以在球体或球壳、圆环之内或之外。

下面分情况讨论。
\subsubsection{均匀球壳与球壳外质点的引力}
因为质点在球壳外，所以r<R,r<=a。均匀球壳密度p=p(r)=常数。引力距离变量u的积分区间为[R-r,r+R]。因为是球形，积分区间为$\phi=[0,\pi],\theta=[0,2\pi]$。于是有

\begin{align}
	F_s&=Gmp\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(\frac{r^2-R^2}{u}+u)|_{R-r}^{r+R}rdr d\theta\\
	F_s&=Gmp\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(r-R+r+R+(r+R)-(R-r))rdr d\theta\\
	F_s&=Gmp\int_0^{2\pi}\frac{1}{2R^2}4r^2dr d\theta\label{MassPointGravityEQ130}
\end{align}

对于球壳，因为密度是面密度p，故厚度可视为单位厚度即dr=1，而球壳半径r=a，代入式\ref{MassPointGravityEQ130}，有

\begin{align}
	F_s&=Gmp\int_0^{2\pi}\frac{1}{2R^2}4a^2 d\theta
	F_s&=Gmp4\pi a^2\frac{1}{R^2}
\end{align}
\subsubsection{均匀球壳与球壳内质点的引力}
因为质点在球壳内，所以a>=r>R。均匀球壳密度p=p(r)=常数。引力距离变量u的积分区间为[r-R,r+R]。因为是球形，积分区间为$\phi=[0,\pi],\theta=[0,2\pi]$。于是有

\begin{align}
	F_s&=Gmp\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(\frac{r^2-R^2}{u}+u)|_{r-R}^{r+R}rdr d\theta\\
	F_s&=Gmp\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(r-R+r+R-(r+R)-(r-R))rdr d\theta=0
\end{align}
\subsubsection{均匀球体与球体外质点的引力}
因为质点在球体外，所以r<R,r<=a。均匀球体密度p=p(r)=常数。引力距离变量u的积分区间为[R-r,r+R]。因为是球形，积分区间为$\phi=[0,\pi],\theta=[0,2\pi]$。于是有

\begin{align}
	F_s&=Gmp\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(\frac{r^2-R^2}{u}+u)|_{R-r}^{r+R}rdr d\theta\\
	F_s&=Gmp\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(r-R+r+R+(r+R)-(R-r))rdr d\theta\\
	F_s&=Gmp\int_0^{2\pi}\frac{1}{2R^2}4r^2dr d\theta\\
	F_s&=Gmp\int_0^{2\pi}d\theta\frac{2}{R^2}\frac{1}{3}r^3|_0^a \\
	F_s&=Gmp\int_0^{2\pi}d\theta\frac{2}{R^2}\frac{1}{3}r^3|_0^a \\
	F_s&=Gmp\frac{4\pi}{3}a^3\frac{1}{R^2}
\end{align}
\subsubsection{均匀球体与球体内质点的引力}
因为质点在球体内，所以a>=r>R。均匀球体密度p=p(r)=常数。引力距离变量u的积分区间为[r-R,r+R]。因为是球形，积分区间为$\phi=[0,\pi],\theta=[0,2\pi]$。于是有

\begin{align}
	F_s&=Gmp\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(\frac{r^2-R^2}{u}+u)|_{r-R}^{r+R}rdr d\theta\\
	F_s&=Gmp\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(r-R+r+R-(r+R)-(r-R))rdr d\theta=0
\end{align}
\subsubsection{均匀薄圆板与薄圆板外质点的引力}
设球半径为a的球体上截取一个薄圆板半径为b，抛弃球体剩余部分，球心O是坐标系原点，在薄圆板圆心C和O的连线上有一个质点A质量为m，OA指向Z轴正方向。OA=R，r<R,r<=a。均匀薄圆板密度p=p(r)=常数。引力距离变量u的积分区间为[R-r,x]。因为是球形，积分区间为$\phi=[0,arcsin(b/a)],\theta=[0,2\pi]$。

现在求u的上限。它是

\begin{align}
	\phi_2&=arcsin\frac{b^2}{a^2}\\
	cos\phi_2&=1-\frac{b^2}{a^2}\label{MassPointGravityEQ132}
\end{align}

代入式\ref{MassPointGravityEQ103}，得

\begin{align}
	u^2&=r^2+R^2-2Rrcos\phi_2\\
	u^2&=r^2+R^2-2Rr(1-\frac{b^2}{a^2})\label{MassPointGravityEQ134}\\
	u^2&=(R-r)^2+2Rr\frac{b^2}{a^2}\label{MassPointGravityEQ136}
\end{align}

于是有

\begin{align}
	F_s&=Gmp\int_0^{2\pi}\int_0^{a}\frac{1}{2R^2}(\frac{r^2-R^2}{u}+u)|_{R-r}^{((R-r)^2+2Rr\frac{b^2}{a^2})^{1/2}}rdr d\theta\\
	F_s&=Gmp\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(\frac{r^2-R^2+(R-r)^2+2Rr\frac{b^2}{a^2}}{\sqrt{(R-r)^2+2Rr\frac{b^2}{a^2}}}+R+r+R-r)rdr d\theta\\
	F_s&=Gmp\int_0^{2\pi}\frac{1}{2R^2}4r^2dr d\theta\\
	F_s&=Gmp\int_0^{2\pi}d\theta\frac{2}{R^2}\frac{1}{3}r^3|_0^a \\
	F_s&=Gmp\int_0^{2\pi}d\theta\frac{2}{R^2}\frac{1}{3}r^3|_0^a \\
	F_s&=Gmp\frac{4\pi}{3}a^3\frac{1}{R^2}
\end{align}
\subsubsection{均匀薄圆环与薄圆环外质点的引力}
因为质点在圆环外，所以r<R,ai<=r<=ao。圆环内径和外径分别是ai,ao,均匀圆环密度p=p(r)=常数。引力距离变量u的积分区间为[R-r,R-r+d]。因为从球形变成圆环，圆环厚度T，则积分区间为$\phi=[0,arccos(2T/(ao+ai))],\theta=[0,2\pi]$。r的积分区间显然是[ai,ao]。于是有

\begin{align}
	F_s&=Gmp\int_0^{2\pi}\int_{ai}^{ao}\frac{1}{2R^2}(\frac{r^2-R^2}{u}+u)|_{R-r}^{R-r+d}rdr d\theta\\
	F_s&=Gmp\int_0^{2\pi}\int_{ai}^{ao}\frac{1}{2R^2}(-(R+r)(1-(d/(R-r))+(R-r+d)+(r+R)-(R-r))rdr d\theta\\
	F_s&=Gmp\int_0^{2\pi}\int_{ai}^{ao}\frac{1}{2R^2}(d(R+r)/(R-r)+1)dr d\theta\label{MassPointGravityEQ130}
\end{align}

\begin{align}
	F_s&=Gmp\int_0^{2\pi}\frac{1}{2R^2}4a^2 d\theta
	F_s&=Gmp4\pi a^2\frac{1}{R^2}
\end{align}
\subsection{均匀球面与球面内质点间的引力}
若质点质量为m，与球心的距离为R。设球的半径为a，R<a，面密度为p，质量为$M=p4\pi a^2$。

直角坐标系和球坐标系坐标原点都在质点，在两个坐标系球心坐标分别为(0,0,R)和(R,0,0)，球面任意微元球坐标为$(r,\phi,\theta)$，建立如图\ref{MassPointGravity2}所示的坐标系。

根据对称性可知，球对质点的引力必沿z方向，x，y方向上合力为0。

微元对质点的引力为dF

\begin{align}
	C&=Gmp\label{MassPointGravityEQ100}\\
	dF&=G\frac{m\rho r^2sin\phi}{u^2}dr d\phi d\theta (R>a)\label{MassPointGravityEQ101}
\end{align}

该引力在z方向上的分力为：

\begin{align}
	dF_z&=dFcos\alpha \label{MassPointGravityEQ102}
\end{align}

由余弦定理，有：

\begin{align}
	u^2&=r^2+R^2-2rRcos\phi\label{MassPointGravityEQ103}\\
	r^2&=u^2+R^2-2uRcos\alpha\label{MassPointGravityEQ104}
\end{align}

式\ref{MassPointGravityEQ103}加式\ref{MassPointGravityEQ104}，得

\begin{align}
	u^2+r^2&=r^2+u^2+2R^2-2rRcos\phi-2uRcos\alpha\\
	0&=R-rcos\phi-ucos\alpha\\
	ucos\alpha&=R-rcos\phi\label{MassPointGravityEQ105}\\
	cos\alpha&=\frac{R-rcos\phi}{u}\label{MassPointGravityEQ106}
\end{align}

对式\ref{MassPointGravityEQ103}两边求微分，得

\begin{align}
	udu&=rRsin\phi d\phi\label{MassPointGravityEQ107}
\end{align}

由式\ref{MassPointGravityEQ103}，得

\begin{align}
	cos\phi&=\frac{r^2+R^2-u^2}{2rR}\label{MassPointGravityEQ108}
\end{align}

式\ref{MassPointGravityEQ106}代入式\ref{MassPointGravityEQ102}，得

\begin{align}
	dF_z&=G\frac{m\rho r^2(R-rcos\phi)sin\phi}{u^3}dr d\phi d\theta\label{MassPointGravityEQ109}\\
	F_s&=\iiint_V dF_z\\
	F_s&=\int_0^{2\pi}\int_0^\pi\int_0^a G\frac{m\rho r^2sin\phi dr}{u^2}\frac{(R-rcos\phi)}{u} d\phi d\theta\label{MassPointGravityEQ110}
\end{align}

式\ref{MassPointGravityEQ107},\ref{MassPointGravityEQ108}代入到\ref{MassPointGravityEQ110}，将自变量$\phi$的函数都替换成u的函数

\begin{align}
	F_s&=Gmp\int_0^{2\pi}\int_0^\pi\int_0^a\frac{r^2sin\phi d\phi dr}{u^2}\frac{(R-rcos\phi)}{u} d\theta\\
	F_s&=Gmp\int_0^{2\pi}\int_0^\pi\int_0^a\frac{rRrsin\phi d\phi dr}{Ru^2}\frac{R-rcos\phi}{u} d\theta\\
	F_s&=Gmp\int_0^{2\pi}\int_0^a\int_{r-R}^{r+R}\frac{rudu dr}{Ru^2}\frac{R-r\frac{r^2+R^2-u^2}{2rR}}{u} d\theta\\
	F_s&=Gmp\int_0^{2\pi}\int_0^a\int_{r-R}^{r+R}\frac{r}{R}\frac{ududr}{u^2}\frac{2R^2-r^2-R^2+u^2}{2Ru} d\theta\\
	F_s&=Gmp\int_0^{2\pi}\int_0^a\int_{r-R}^{r+R}\frac{r}{2R^2}\frac{R^2-r^2+u^2}{u^2}dudr d\theta\\
	F_s&=Gmp\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(\frac{r^2-R^2}{u}+u)|_{r-R}^{r+R}rdudr d\theta\\
	F_s&=Gmp\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(r-R+r+R-(r+R)-(r-R))rdudr d\theta=0
\end{align}
\subsection{均匀球面与球面外质点间的引力}
设质点质量为m，与球心的距离为R。球的半径为a，R>a，球面密度为$\rho$，质量为$M=\rho4\pi a^2$。直角坐标系和球坐标系坐标原点都在球心，在两个坐标系质点坐标分别为(0,0,R)和$(R,0,0)$，球内任意微元球坐标为$(r,\phi,\theta)$，其体积为$dV=r^2sin\phi dr d\phi d\theta$，微元和质点的距离为u，建立如图\ref{MassPointGravity}所示的坐标系。

\begin{figure}
	\centering
	\epsfig{file=MassPointGravity,width=8cm}
	\caption{质点与均匀球面间的引力\label{MassPointGravity}}
\end{figure}

根据对称性可知，球对质点的引力必沿z方向，x，y方向上合力为0。

微元对质点的引力为dF

\begin{align}
	C&=Gm\rho\label{MassPointGravityEQ100}\\
	dF&=G\frac{m\rho r^2sin\phi}{u^2}dr d\phi d\theta \label{MassPointGravityEQ101}
\end{align}

该引力在z方向上的分力为：

\begin{align}
	dF_z&=dFcos\alpha \label{MassPointGravityEQ102}
\end{align}

由余弦定理，有：

\begin{align}
	u^2&=r^2+R^2-2rRcos\phi\label{MassPointGravityEQ103}\\
	r^2&=u^2+R^2-2uRcos\alpha\label{MassPointGravityEQ104}
\end{align}

式\ref{MassPointGravityEQ103}加式\ref{MassPointGravityEQ104}，得

\begin{align}
	u^2+r^2&=r^2+u^2+2R^2-2rRcos\phi-2uRcos\alpha\\
	0&=R-rcos\phi-ucos\alpha\\
	ucos\alpha&=R-rcos\phi\label{MassPointGravityEQ105}\\
	cos\alpha&=\frac{R-rcos\phi}{u}\label{MassPointGravityEQ106}
\end{align}

对式\ref{MassPointGravityEQ103}两边求微分，得

\begin{align}
	udu&=rRsin\phi d\phi\label{MassPointGravityEQ107}
\end{align}

由式\ref{MassPointGravityEQ103}，得

\begin{align}
	cos\phi&=\frac{r^2+R^2-u^2}{2rR}\label{MassPointGravityEQ108}
\end{align}

式\ref{MassPointGravityEQ106}代入式\ref{MassPointGravityEQ102}，得

\begin{align}
	dF_z&=G\frac{m\rho r^2(R-rcos\phi)sin\phi}{u^3}dr d\phi d\theta\label{MassPointGravityEQ109}\\
	F_s&=\iiint_V dF_z\\
	F_s&=\int_0^{2\pi}\int_0^\pi\int_0^a G\frac{m\rho r^2sin\phi dr}{u^2}\frac{(R-rcos\phi)}{u} d\phi d\theta\label{MassPointGravityEQ110}
\end{align}

式\ref{MassPointGravityEQ107},\ref{MassPointGravityEQ108}代入到\ref{MassPointGravityEQ110}，将自变量$\phi$的函数都替换成u的函数

\begin{align}
	F_s&=Gm\rho\int_0^{2\pi}\int_0^\pi\int_0^a\frac{r^2sin\phi d\phi dr}{u^2}\frac{(R-rcos\phi)}{u} d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^\pi\int_0^a\frac{rRrsin\phi d\phi dr}{Ru^2}\frac{R-rcos\phi}{u} d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\int_{R-r}^{r+R}\frac{rudu dr}{Ru^2}\frac{R-r\frac{r^2+R^2-u^2}{2rR}}{u} d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\int_{R-r}^{r+R}\frac{r}{R}\frac{ududr}{u^2}\frac{2R^2-r^2-R^2+u^2}{2Ru} d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\int_{R-r}^{r+R}\frac{r}{2R^2}\frac{R^2-r^2+u^2}{u^2}dudr d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(\frac{r^2-R^2}{u}+u)|_{R-r}^{r+R}rdr d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(r-R+r+R+(r+R)-(R-r))rdr d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\frac{2}{R^2}a^2d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}d\theta\frac{2}{R^2}a^2\\
	F_s&=Gm\rho\frac{4\pi a^2}{R^2}
\end{align}

所以均匀球体对球外一点的引力好象球体的质量全部集中在球心一样。那么两个均匀球体间的引力就可以分别把质量全部集中至各自球心，所以用公式计算时R就是球心间距离。
\subsection{均匀薄圆环与环内质点的引力}
\begin{figure}
	\caption{圆形薄片引力\label{circlethinpinGravity}}
	\includegraphics[scale=0.3]{circlethinpinGravity}
\end{figure}
设质点质量为m，与环心的距离为R。球的半径为a，R>a，球面密度为$\rho$，质量为$M=\rho4\pi a^2$。直角坐标系和球坐标系坐标原点都在球心，在两个坐标系质点坐标分别为(0,0,R)和$(R,0,0)$，球内任意微元球坐标为$(r,\phi,\theta)$，其体积为$dV=r^2sin\phi dr d\phi d\theta$，微元和质点的距离为u，建立如图\ref{MassPointGravity}所示的坐标系。

\begin{figure}
	\centering
	\epsfig{file=MassPointGravity,width=8cm}
	\caption{质点与均匀球面间的引力\label{MassPointGravity}}
\end{figure}

根据对称性可知，球对质点的引力必沿z方向，x，y方向上合力为0。

微元对质点的引力为dF

\begin{align}
	C&=Gm\rho\label{MassPointGravityEQ100}\\
	dF&=G\frac{m\rho r^2sin\phi}{u^2}dr d\phi d\theta \label{MassPointGravityEQ101}
\end{align}

该引力在z方向上的分力为：

\begin{align}
	dF_z&=dFcos\alpha \label{MassPointGravityEQ102}
\end{align}

由余弦定理，有：

\begin{align}
	u^2&=r^2+R^2-2rRcos\phi\label{MassPointGravityEQ103}\\
	r^2&=u^2+R^2-2uRcos\alpha\label{MassPointGravityEQ104}
\end{align}

式\ref{MassPointGravityEQ103}加式\ref{MassPointGravityEQ104}，得

\begin{align}
	u^2+r^2&=r^2+u^2+2R^2-2rRcos\phi-2uRcos\alpha\\
	0&=R-rcos\phi-ucos\alpha\\
	ucos\alpha&=R-rcos\phi\label{MassPointGravityEQ105}\\
	cos\alpha&=\frac{R-rcos\phi}{u}\label{MassPointGravityEQ106}
\end{align}

对式\ref{MassPointGravityEQ103}两边求微分，得

\begin{align}
	udu&=rRsin\phi d\phi\label{MassPointGravityEQ107}
\end{align}

由式\ref{MassPointGravityEQ103}，得

\begin{align}
	cos\phi&=\frac{r^2+R^2-u^2}{2rR}\label{MassPointGravityEQ108}
\end{align}

式\ref{MassPointGravityEQ106}代入式\ref{MassPointGravityEQ102}，得

\begin{align}
	dF_z&=G\frac{m\rho r^2(R-rcos\phi)sin\phi}{u^3}dr d\phi d\theta\label{MassPointGravityEQ109}\\
	F_s&=\iiint_V dF_z\\
	F_s&=\int_0^{2\pi}\int_0^\pi\int_0^a G\frac{m\rho r^2sin\phi dr}{u^2}\frac{(R-rcos\phi)}{u} d\phi d\theta\label{MassPointGravityEQ110}
\end{align}

式\ref{MassPointGravityEQ107},\ref{MassPointGravityEQ108}代入到\ref{MassPointGravityEQ110}，将自变量$\phi$的函数都替换成u的函数

\begin{align}
	F_s&=Gm\rho\int_0^{2\pi}\int_0^\pi\int_0^a\frac{r^2sin\phi d\phi dr}{u^2}\frac{(R-rcos\phi)}{u} d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^\pi\int_0^a\frac{rRrsin\phi d\phi dr}{Ru^2}\frac{R-rcos\phi}{u} d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\int_{R-r}^{r+R}\frac{rudu dr}{Ru^2}\frac{R-r\frac{r^2+R^2-u^2}{2rR}}{u} d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\int_{R-r}^{r+R}\frac{r}{R}\frac{ududr}{u^2}\frac{2R^2-r^2-R^2+u^2}{2Ru} d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\int_{R-r}^{r+R}\frac{r}{2R^2}\frac{R^2-r^2+u^2}{u^2}dudr d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(\frac{r^2-R^2}{u}+u)|_{R-r}^{r+R}rdr d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\int_0^a\frac{1}{2R^2}(r-R+r+R+(r+R)-(R-r))rdr d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}\frac{2}{R^2}a^2d\theta\\
	F_s&=Gm\rho\int_0^{2\pi}d\theta\frac{2}{R^2}a^2\\
	F_s&=Gm\rho\frac{4\pi a^2}{R^2}
\end{align}
